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Q I We consider a simple model of the dynamics of a single electron 

Q ■ in a crystal lattice. Although this is a standard problem in condensed 

. matter physics, alternative ways of evaluating a partition function for 

such a system lead to equalities, that may be interesting from the point 
(—1 ' of view of mathematical analysis, combinatorics and graph theory. 
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In this paper we present results which may be interesting from the point of 
view of pure mathematics, but originate from some properties of mathemat- 
ical formalism used to describe certain physical systems. Namely, we derive 
some generalizations of the integral representation of the modified Bessel 
^ . function Jq (2^), 

■ ^ I „2£cos(k)jK ST^ ^ 
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here given in the form of the Taylor expansion around ^ = (see pD], 9.6.16, 
^ \ p. 376). The important observation is that formula ([I]) describes a 'partition 

\ junction of the tight-binding model of a single electron moving in an infi- 

nite, one-dimensional chain of atoms. This model and its generalizations are 
widely used in condensed matter physics and may be found in any standard 
textbook on that field, (see e.g. [2], p. 11). 

If we consider, in the framework of the same model, some D-dimensional 
crystal lattice, and if we write down the corresponding partition function, we 
obtain equality analogous to ([T]). Namely, its l.h.s is then a sum or integral 
over some subset of M^, having nontrivial symmetry properties. On the other 
hand, the r.h.s. of both formula ([T]) and any of its generalizations is a series 
in one or several complex variables, and can be found with the help of the 
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method described in detail in the present paper. In that way, we are able to 
establish several identities, which may be interesting from the point of view of 
mathematical analysis. What is more important, our method involves some 
basic techniques and concepts of combinatorics and graph theory, therefore 
it may be interesting from the point of view those disciplines. 

The results presented here have been obtained with the help of a standard 
mathematical formalism and techniques used in theoretical condensed-matter 
physics. This fact has two important consequences. First, our notation is 
generally an usual one for that field. This results in occurrence of additional 
numerical factors and indices, which may be found redundant, in the formu- 
las presented. Hence, throughout the paper we use units in which Planck 
constant h and Boltzmann constant ks are equal unity, h = ks = 1. 

Second, our results may be found not rigorous enough from the point of 
view of pure mathematics. Nonetheless, we believe that they are correct, and 
they can be given a more formal shape if necessary. This should be one of 
the subjects of our future investigation. 

This paper is organized as follows: in section 1. we introduce some con- 
cepts necessary to understand our method and we describe the method itself. 
In section 2. we analyze various lattice geometries and present the results 
obtained. Section 3. consists of appendices, where some details, omitted in 
the main text, are presented. 

1 The method 

1.1 Preliminaries: physical context 

Let us consider a motion of a single electron in a periodic crystal lattice. To 
construct mathematical model of such a physical system, we have to provide 
a description of a crystal lattice, as well as to express the quantum dynamics 
of the electron on that lattice. 

To model a crystal lattice, we define first the D-dimensional abstract 
lattice, A, in a standard manner, as a set of lattice vectors R, r G M'^, 
being the linear combinations with integral coefficients of D fundamental 
translation vectors ai, a2, . . . , a/) 



However, appropriate model of the real crystal structure is constructed by 
attaching basis of M > 1 atoms to each point of the abstract lattice (see 
e.g. [3], p. 5). We denote positions of atoms by vectors i, called lattice sites, 
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each lattice site has its corresponding lattice vector, R(i). For lattices with 
M = 1, lattice vectors can be chosen in such a way that they coincide with 
lattice sites. This is the case of Bravais lattice, and most, although not all, 
cases analyzed in this paper belong to this class. On the other hand non- 
Bravais lattices consist of M > 1 pairwise disjoint sublattices, with any two 
points of a given sublattice being equivalent. 

It is of crucial importance not to confuse A with the set of lattice sites 
i, denoted H. From now on we use the term lattice only for S, whereas 
A is always termed abstract lattice, as before. A = |A| (S = |H|) is the 
number of the lattice vectors (sites), respectively. Obviously, S = MA, and 
this number may be finite or not. I the former case the periodic boundary 
conditions (PBC) have to be imposed to ensure periodicity. 

Apart from abstract lattice A we have to remind the Reader also the 
notion of the reciprocal lattice, defined as the set of vectors K having the 
following property (see e.g. [3], p. 61), 

giK.R _ R G A. (3) 

Now we return to the question of a quantum- mechanical description of 
the dynamics of the electron. To be able to give such a description, first we 
have to specify the Hilbert space 7i of the state vectors, and then choose some 
basis vector set in this space. It turns out that set H in a real space provides 
the convenient basis in Ti. Namely, we define state |i) G 7i as describing 
the electron cloud centered on the lattice site i. This is so called Wannier 
basis or Wannier representation 0. Such defined states form a complete and 
orthonormal basis set, 

5^|i)(i|=i, (i|j) =S,y (4) 

iGS 

It is convenient to treat our Hilbert space as a subspace of the larger one, 
called the Fock space, T . This new space, apart from the just defined one- 
particle states |i) , consist of other iV-particle states, with iV = 0, 2, . . . , Nmax, 
Nmax = and may be written as a direct sum of its A^-particle sectors, = 
J^o © -^1 © ■ ■ ■ © ^N^ax ■ Although we are interested in one-particle sector only, 
J^i = Ti, we use the general formalism of second quantization, as the most 
convenient for our purposes (compare chapter 1). We introduce electron 
creation (annihilation) operators c\ (cj) for lattice site i (j), respectively. 

^It may be also regarded as some kind of coordinate representation of quantum me- 
chanics, but in the discretized, rather then continuous, real space. We use those three 
terms interchangeably in the following sections. 
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They obey standard anti-commutation relations 

c-A + c\cj = 6ij, CjCi + CiCj = 0, (5) 

and are related to the basis state vectors |i) ((i|) in the following way 

|i) = cl|0), (i| = (0|ci, (6) 

where |0) is a vacuum state, spanning zero-particle sector JFq C J-'. 

The dynamics of the electron can be then described as follows: initially 
localized in the vicinity of a particular site i, it may quantum-mechanically 
tunnel to neighboring site j. In 7i it corresponds simply to the transition 
from basis state |i) to state |j). If we include only the so-called single-band 
dynamics and neglect other possible degrees of freedom (like spin, orbital, 
etc.) of the electron, the resulting Hamilton operator (Hamiltonian) reads 

H=J2tih})4ci= Yl ^(q)4ci+q- (7) 

iJes ieH,qe5 

Here the summation over i runs through the whole lattice H. We assume that 
the complex numbers t(i, j) {hopping integrals) fulfill the following conditions, 

Vijes: t(i,j)=t(j,i)* = t(j-i)=t(q), q = j - i, (8) 

due to the fact that Hamiltonian is a hermitian operator, and because all 
lattice sites of a given sub-lattice are assumed to be equivalent [translational 
invariance) . We also put t{i, i) = for all i. By S we denote a set of vectors q, 
for which t(q) ^ 0. We assume the maximal number of independent hopping 
integrals t(q) = € {^1,^2, • • • ,tc}, namely C = the factor ^ appears 

due to (IH])- Consequently, in general, the Hamiltonian does not possess all 
symmetry properties of the lattice - the discrete rotational symmetry may 
be lost unless we impose specific conditions on parameters tg. 

Due to translational invariance of the Hamiltonian ([7]), it is advisable to 
express it in Block basis, labeled by the values of quasimomentum k, and 
related to the Wannier basis by the unitary transformation, 

cl\0) ^ |k) = ^ 5^ exp ( - zk ■ R(i)) |i) (9) 

with the inverse given by 

cj|0) = |i) = ^ E exp(zk-R(i))|k). (10) 



4 



The above formulas refer only to the case of finite (A < oo), Bravais lattice 
(M = 1). The formalism for non-Bravais lattices is not presented here, the 
limit A — s> oo will be considered later. 

The vectors k in ( PfTOl) . determined by the lattice geometry and PBC 
imposed, belong to the first Brillouin zone (FBZ), which is the Wigner-Seitz 
unit cell of the reciprocal lattice centered at K = (compare [3J, chapter 9). 

Formulas ^ and (ITUI) give us also the transformation properties 
of the creation (annihilation) operators themselves, and can easily find the 
form of the Hamiltonian (JTj) in the Bloch basi^, in which it is diagonal, 
-ff|k) = e(k)|k). Explicitly, we have 

^ = ^e(k)4ck = J]e(k)nk. (H) 

k k 

The functional dependence of its eigenvalues e(k) on k, the dispersion rela- 
tion, results from the lattice geometry, but is also a function of the hopping 
integrals tg-, e = e(k; ti, . . . , tc)- For M = 1 it has the following general form 

e(k; ti, . . . , tc) = tiei(k) + . . . + tcec(k). (12) 

We assume, that e(k; ti, . . . , tc) is at least a piecewise continuous, and then 
integrable, function of k. 

Consider now a limit A — >■ cxd. Due to the periodic boundary conditions 
(PBC), for any finite A the number of k vectors in FBZ is equal that of 
lattice vectors R, R G A. As a consequence, in the limit in question, some 
physical quantities, expressed as sums over k G FBZ, may be divergent. To 
keep the values of those quantities finite, we rescale them by the factor A, 

E /(k) — i E /(k)- (13) 

The sum on the r.h.s of ( fT3l) . for a given value of A, may be regarded as 
the Riemann sum of a function /(k). The important fact is that a mesh 
of any partition of FBZ given by the values of k goes to zero in the limit 
A — > CX3, regardless of the particular choice of PBC. Consequently, in that 
limit, (fT3l) approaches the value of the Riemann integral over a FBZ of /(k), 
here understood as a function of continuous variable k. We can then replace 
the summation in ffT^ by the integration over FBZ, according to 

where ^ is a Jordan measure of the first Brillouin zone. 



^Referred to also as the quasimomentum representation or h- representation. 
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1.2 The partition function 



1.2.1 General formalism and quasimomentum representation. 

According to quantum statistical mechanics (see e.g. [5j, p. 245), the parti- 
tion function, is given by the formula 

Z = Z(/?) = Tr(exp(-/?^)), (15) 

where Tr denotes trace of an operator, (3 = 1/T is an inverse temperature 
and H is the Hamiltonian of a system considered. We evaluate Z first in 
Bloch basis, which is rather simple, due to the fact that one-electron Bloch 
states cj^|0) are eigenstates of the Hamiltonian (JTj). In that case the partition 
function reads 

FBZ 

for finite number of lattice vectors, and 

Z = i / exp(-/5e(k))cik (17) 

^ J FBZ 

for the A — i> oo limit. To avoid infinite value of Z in ( fT6i) we have applied nor- 
malization described in a previous paragraph by formula (fT3l) . Introducing 
for convenience C new variables C,s, defined as follows 

is = -f3t,, sG{l,2,...,C}, (18) 

and having in mind (fT2!) . we can rewrite both (fT6l) and (ITTl) as 

1 ^ 

Z = ^(6,6, ■■■,^c) = jYl exp (5^ese.(k)) (19) 

FBZ 1 

and 

(J 

Z = Z(6,6,...,ec) = :^ / exp(Ve.e,(k))dk. (20) 

^ J FBZ ^ 

Let us point out, that for most lattice geometries considered, the sum 
(fT9ll or the integral (l20l) cannot be easily evaluated, and the explicit form 
of the functional dependence of Z{S,i, ^2, • • • , ic) on its arguments cannot be 
evaluated either in terms of elementary functions or in terms of standard 
special functions. 
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1.2.2 Coordinate representation 

The general formula for the partition function of our one-electron system, 
evaluated with the help of Bloch basis states ( |T9l) . (l20l) . has a form of a 
complicated sum or integral, which usually cannot be handled analytically. 
However, the trace of an operator does not depend on the basis one uses, and 
in order to obtain an alternative expression for Z{^i, ^2, • • • , ^c) "we compute 
it now using Wannier basis. 

Evaluation of the partition function in the present case requires some 
effort, below we describe a method used for this purpose. First, assuming 
A < 00, we rescale Z by A, in correspondence with the formula f|T3|) . 

iec ies n 

The problem in question is solved, once we know how to evaluat^ the quan- 
tity g{n) defined as 

l5^(i|(-/?^r|i) ^g{n), (22) 
ies 

for the Hamiltonian (171). Explicitly, in this particular case we have 

9{n) = lj2(^\{ (23) 
ies i,j(i)es 

where we put t(ij) = and introduce a shorthand = cjcj to make the 
following formulas more compact and transparent. One can easily examine 
some obvious properties of just introduced operators Aij, below we invoke 
only two of them, 

A\. = Ay„ AjAji = ni(l-nj). (24) 
When calculating g{n) we have to deal with the terms of the general form 

(~/^)"^l2nl2n-l'^l2n-2l2n-3 ' ' ' '^Uls ^bU ^l2nl2n- 1 ^l2n-2l2n-3 ' ' ' ^Uh^hhj (^5) 

where Im G S. In what follows we concentrate on the operator part of ( l25l) . 

^l2nl2n-l^l2n-2l2n-3 " " " -4l4l3 ^l^l^ . (26) 

^We postpone for a moment a question if it is legitimate to interchange the summation 
over the order of the expansion n and that over the states |i). 
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Most of such terms do not contribute to g{n). This is, firstly, due to the 
fermion anticommutation relations and properties of the occupation number 
operatoi0, that all sequences of the type 

■■■AiiAij--- , ■■■AiiAji--- or ■■■ni(l-ni)--- , (27) 

identically vanish. Secondly, because we compute trace in basis of one- 
electron states, all terms where the product of two or more occupation elec- 
tron number operators for different sites, i, j appear, 

ni(l - ni)nj(l - nm) ■ ■ • , (28) 
must vanish, too. This restricts the form of the terms in question to 

An+lin^nin-l " " " ^1312^211 = " ^in ) ( 1 " ^in-l) " " " (1 " ^12)^1- (29) 

However, due to the trace operation we have to consider only diagonal terms, 
i.e. those with the first and the last index in any product of n Ai^ operators 
being identical. Consequently, the only allowed terms in an expansion giving 
nonzero contribution to g{n) must be of the following form 

AlinAnin-l " ■■^312^211 = ^1 ( ^ " ^in)(l " ^in-l) " " " (1 " ^iJCi,. (30) 

Notice that the above expression is a product of n — 1 hole number operators 
sandwiched in between c| and operators, and has trace equal unity. 
Now let us now return to the full form of the terms ( l30i) . 

(— /?) i!^iiin^inin-l '^1312^1211 ^ilin^inin-l ^1312^1211 ) (^l) 

and assume for a moment that among n hopping integrals in ( 13T1) , tg appears 
Hs times, s = 1, 2, . . . , C. Having in mind variables introduced in (ITSjl . we 
see, that (IHTj) contributes to g{n) given by (1^ through the term 

SI S2 ■ ■ ■ SC ' 

We denote the total number of terms contributing in the same way to g{n) 
by r(ni, n2, . . . , nc)- To obtain g{n) we must sum up the terms (132|) for all 
permitted values ni, n2, . . . , ric, 

constr 

g{n; ei, 6, • • • , ec) = Yl ^2, . . . , nc) ■ er^s^ ■■■Cc^, (33) 

ni,n2,...,nc 

^More precisely, we refer to fii as to an electron and 1 — fii as to a hole number operator 
for a given site i, respectively. 



C 



n. 



(32) 
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where we have exphcitly denoted the dependence of g{n) on .^1,^25 • • • ,C,c, 
and the fact that summation in flH^ obeys the constraint appearing in 
Comparing (EH), (El and ([33]), we have 



constr -p/ N 

^(ei,6,...,ec) = E E (34) 



n ni,n2,...,nc 



Obviously, from the above formula it follows that the r(ni, n2, . . . , ric) coef- 
ficients are related to the respective coefficients Zn-^n2...nc of the Taylor series 
expansion of ^2, • • • , ^c), 

r(ni, 77,2, ...,nc)=n\- Zn^n2...nc, (35) 

defined as 

gni Qn2 Que Z{i^, . . . 



Zn,n2...nc Q^n, Q^n2 QQc ^^j^^j . . . ^^j 



■ (36) 

5l=6=-=«C=0 



However, according to the remark in section 2.1, in most cases fll9p . (1201) ) and 
their derivatives with respect to cannot be evaluated analytically. In order 
to find an explicit formula for the Taylor series of the partition function, we 
have to develop an alternative method of computing coefficients Zn^ri2...nc (o^' 
r(?7,i, n2, . . . , nc))- This problem is solved in the next chapter. 



1.2.3 The method of evaluating r(rai,n2, . . . ,nc)- 

There exist a method of combinatorial character, which allows us to find 
r(ni,n2, . . . ,nc) coefficients. It is based on the following observation: the 
dynamics of our system, governed by the Hamiltonian (I7j) may be described 
by an intuitive picture of an electron tunneling - or 'taking a step' form one 
site to the other. More precisely, such a step is generated by Aij operator, and 
may be visualized by an arrow, pointing from site j to site i. Each such arrow 
is labeled by a complex number - a hopping integral ty = tg. Consequently, 
the product of operators given by (IHUI) can be visualized, using terms of 
graph theory, as a closed path or a cycle^ with the respective lattice sites 
being its vertices, and with the complex number ts ascribed to each edge. 
The crucial fact is the one-to-one correspondence between terms fl30l) and 
the resulting closed paths, for details see appendix 13. II For our purposes it is 
important which vertex of our graph is a terminal one, (i.e. this from which 
our path starts and where it ends), the graphs differing by a terminal vertex 
are different from our point of view. Notice also, that any vertex, including a 
terminal one, may be 'visited' by an electron more then once, i.e. there may 
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be more than one edge connecting two vertices, and thus our closed paths 
does not have to be simple. 

The set of all vertices that can be connected by a directed edge with a 
particular initial site are given by vectors q G 5. In other words, electron 
can move in direction and at the distances determined by those vectors, or 
equivalently, by the hopping integrals t(q). (This is the reason why it is 
convenient to assume the maximal number of distinct hopping integrals, i.e. 
one t for every pair (q, — q), ±q G S). The above described crystal directions 
should not be confused with D mutually orthogonal spacial directions, which, 
in principle can be chosen at will. Obviously, for the closed path we have 
to fulfill the condition that a D-dimensional vector of the total displacement 
must be equal zero. 

At first glance it seems than in order to find r(?T,i, n2, . . . , nc) for a given 
lattice geometry, we have to compute the total number of the closed paths of 
lengtho n, with exactly Ug edges labeled by the respective hopping integral ts- 
However, due to the translational invariance of the Hamiltonian and our nor- 
malization procedure, the summation over terminal sites in trace operation 
just cancel the normalization factor 1/A. What remains then is to compute 
the number of such paths with one particular site chosen as a terminal one. 

Number of closed paths of length n, with Ug edges 
r(ni, n2, . . . , uc) = labeled by tg, and with one particular site i 

chosen as the terminal vertex. (37) 

This allows us to compute Z(^i, • • • , ^c) also in the A ^ oo limit, the 
difference with the case of finite lattice is that in the latter, when computing 
r(ni, n2, . . . , nc), periodic boundary conditions must be taken into account. 

The number of independent variables in Z{^i, ^2, • • • , ic) can be reduced 
from its maximal value C by putting = for some s, w, after all coef- 
ficients r(ni,n2, . . . ,nc) are found. In the extreme case, partition function 
becomes function of one variable ^ only, and the total number of all closed 
paths of the total length n, starting from a particular lattice site, denoted 
by T{n), is equal 

constr 

f (^) = r(ni, n2,..., nc) =^(n;l, !,...,!). (38) 

ni,n2,...,nc 

To summarize the above discussion, we make the following observation. 
Define a graph jiH), such that its vertices are lattice sites of the lattice 

^By length we mean, according to the terminology of graph theory, the number of edges, 
i.e. the order of the expansion of exp(— and not the distance in a real space. 
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considered, and that vertices j and i are connected by an edge if and only 
if tij 7^ 0. In the case of finite lattice and for all hopping integrals equal 
unity, tij = t = 1, the matrix of the Hamiltonian in the Wannier represen- 
tation is simply the adjacency matrix of the graph 7(i^). This observation 
immediately shows the essence of our method, the n-th power of adjacency 
matrix contains information about the number of all paths of length n, its 
trace gives the total number of all closed paths. Generalization to the case 
of complex tg requires some kind of weighted graphs. 

For infinite lattices such a simple picture is no longer available. 

1.2.4 The main formula 

I the previous section we have presented the method of finding the coef- 
ficients r(ni,n2, . . . ,nc), which does not require exphcit evaluation of the 
sum ( HM or integral (pUj) nor the respective partial derivatives. Equating 
two alternative expressions for Z, formulas (fT9l) or (!20l) and (jMl) . we obtain 

-, C constr t-\ f / \ 

FBZ 1 n ni,n2,...,nc 

(39) 

or, for the case of infinite lattic^, 

-. » C constr -pj/ x 

i/ exp(^6e,(k))dk = ^ J: ''"""^.-""W -g. 

J FBZ ^ ni,n2,...,nc 

(40) 

This is the main result of the present paper. Above formulas, for particular 
cases of lattice geometries, generate mathematical identities that may be 
even yet unknown. 

There are two points of view one can have on flH^ . fj^Uj) . Firstly, its r.h.s is 
an explicit expression for the series expansion of Z[^i, ^2, ■ ■ ■ , ^c) given by the 
nontrivial sum or integral appearing on the l.h.s. Thus, it could be regarded 
as some kind of generalization of the well-known formula ([T]) giving a modified 
Bessel function Jo(2^). Although the occurring series are rather complicated 
and may be found inconvenient for practical computational purposes, they 
may be interesting from the point of view of the mathematical analysis. 

On the other hand, quasimomentum representation of Z {(,1,(^2, ■■■ ,^c) 
may be used to obtain a coefficients r(ni,ra2, . . . ,nc) - numbers of closed 

^Throughout the text we keep the same symbol F for both A < 00 and A — > 00, only 
here we distinguish both cases by respective superscript to avoid misunderstanding. 
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paths of some specific kind, for a given lattice geometry. Namely, let us 
remind the Reader, that according to and (HUj) . flHUl) we have 



n(ni,n2, ...,nc) = — r-; — r a^n^ a^«2 J^THT , (41) 



n! '9"(^/^^^exp(^f e.es(k))rfk^ 



r(ni,n2,...,nc) = — j — r ;icnif>cn2 ' (^2) 

for finite and infinite lattice, respectively. The r.h.s of the above formulas can 
be easily computed numerically, and therefore our results may be interesting 
from the point of view of graph theory or combinatorics. 

Finally, in the case of all equal, due to (1381) formula (HOl) reads 



(43) 



m = U exp(^5:e.(k))rfk = 5:Me 

V Jfbz 1 „ ri\ 

In the next section the above equality is utilized frequently instead of the 
most general case 



2 The results 

In this part we present our results. For each case we begin with the brief 
description of the geometry of the direct as well as the reciprocal lattice. 
Then, the set iS, as well as C- the number of independent hopping integrals 
in the Hamiltonian (I7j), the dispersion relation and the partition function in 
quasimomentum representation are given. Next, we explain in detail how to 
compute the partition function in Wannier representation, for each partic- 
ular lattice geometry. Finally, we compute explicitly first few terms of the 
Taylor series expansion of Z, both by combinatorial techniques as well as by 
numerical integration of the partition function in k-representation. 

Except for the linear chain lattice (see l2.1.1|) . we consider infinite lattices 
only, and most analyzed cases are Bravais lattices, except for paragraph 12.41 
where we deal with the non-Bravais honeycomb and diamond lattices. Real 
hopping integrals are assumed throughout the main text, the case of complex 
integral for particular lattice geometry is analyzed in appendix 13.21 
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2.1 Linear chain with nearest neighbor 
hopping 

In the present paragraph we consider the hnear chain of atoms, with the 
Hamihonian ([7]) given by 

i q=-l,l 

i.e. with the only nonzero tg between nearest neighboring sites. 
2.1.1 The case of a finite lattice 

We begin from the case of finite number of sites, consequently, we have 
to impose the periodic boundary conditions (PBC). Diagonalization of the 
Hamiltonian (jHj) by Fourier transform yields a dispersion relation 



efc = 2t cos (A;) , (45) 
where due to the PBC quasimomentum k is given by 

k = m e Z. (46) 

To find A allowed values of m (or, equivalently, k), i.e. the first Brillouin 
zone (FBZ), we consider separately the case of even and odd A. We obtain, 
respectively 

mG + A even, (47) 

, (A-1) (A-1) (A-1), , 

m G {--^, + 1, • • • , ^^}' ^ odd. (48) 

The resulting values of k are then 

, , 7r(A-2) 7r(A-2) , 

k G { . . . , tt}, a even (49) 

, 7r(A-l) 7r(A-3) 7r(A-l), . , 
k G {-^^, • • • , ^ ^ A odd. (50) 

Obviously, the permitted values of k are related to the A roots of equation 

= 1, zeC (51) 

in the following way 

k = axgz. (52) 
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The partition function depends now on one variable, = Zc{^), ^ = —fit, 
expressed in quasimomentum representation (fTTl) it reads 

zae,A) = ^ E (53) 

keFBZ 

with the FBZ given by (H9 l 1501) . We denoted exphcitly the dependence of Zc 
on both its natural argument ^ as well as on the number of lattice sites A. 

In order to evaluate Zc{C,,A.) in the Wannier representation, notice, that 
from any site electron can move to the nearest neighboring site either in 
positive (clockwise) or negative (anti clockwise) direction along the chain. 
Denote number of such steps by and n~ , respectively. Obviously, n = 

+ n~ , we also define 

d = d{n,A) = - . (54) 

For a closed path of lenght n on the lattice with A sites, the following con- 
dition has to be fuUfiled 

d{n,A) = cA, ceZ. (55) 

In order to find the permitted values of c, let us note, that for any n and A 
we have 

n = wA + r, r < A. (56) 

or, equivalently, 

U7A = n - n(modA) = A[n/AJ. (57) 

Formula (|56|) shows, that the electron can wind at most w times around the 
chain, and then have still r steps to make to complete the path. Consequently, 
c can take any value between —w and w, 

c G {—w, —w + 1, . . . ,w — l,w} = Cc, (58) 

but r must be an even number, and for even A only the even powers of ^ are 
present in the expansion of Zc{^). 

The number of closed paths of lenght n with exactly steps in the 
positive direction is (J^) = n\/{n~^\n~\), the summation over all allowed 
values of yields their total. 



r,(n;A)_X;^^,^_, + (59) 
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due to fl3^ and where values of d (or n+) are given by fl3^ [3H|) . Finally, the 
partition function reads 



n=0 n=0 d(n) 



n=0 



Equating ( l60l) and ( l53i) . having in mind ( H6l) . we obtain 

^^(^'^) = A E^^P(2^^o<^)) =E^"(^)^"' (61) 

m(A) n=0 

the values of m(A) are given by (H71) and (H8l) . 



2.1.2 The case of infinite lattice 

Consider now the limit A ^ oo. Then no path can be closed by moving 
around the chain. Consequently, only terms with c? = contribute to fl3IH) . 
( l60l) and n must be even. Defining: v = n/2, we can write ( l60l) as 

n=0,even ^ ' ^ ^ ' ^ i/=0 

In the limit considered FBZ = [— vr, vr] and we replace the summation over 
FBZ by an integration (compare formula ( |T^ ). Equating two formulas for 
Zc{C,), analogously to (!6T|) we find 

This is the equality ^ quoted in the introduction, and Zc{i) = /o(20- 

2.2 Linear chain with nearest and next-nearest neigh- 
bor hopping 

Also in this paragraph we analyze the case of an infinite linear chain, so 
the direct and reciprocal lattice geometry are the same as in section 3.2. 
However, now we assume slightly more general form of the Hamiltonian, with 
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the nonzero values of the nearest- as well as next-nearest-neighbor hopping 
integrals, 

Hl = ^Y^ ihclci+, + t2c\ci+2q). (64) 

i q=-l,l 

This difference results in a change in dispersion relation, 

eL{k) = 2ti cos{k) + 2^2 cos(2A;), k E [-vr, tt] = FBZ. (65) 

In the present case the partition function depends on two variables = — 
^2 = — /3t2- Written in the quasimomentum representation it reads 

^l(6,6) = ^ f e2«i™^W+2«^™^(2fe)^^ ^ ^ (66) 

Let us now turn to the problem of evaluating ^^(,^1,^2) in coordinate repre- 
sentation. We start with the simple observation that there are two types of 
'steps', related to hopping integrals ti and ^2- Their length in the real space, 
measured in lattice constants, is equal 1 and 2, respectivelj|l|. Let ni and ^2 
be the numbers of corresponding steps, with n = rii + n2 being their total. 
In contrary to the situation in paragraph 3.2, in the present case n can be 
odd, but rii cannot, as it is impossible to close a path in such Con- 
sequently, n and n2 must have the same parity, however for rii = 0, n2 = n 
must be even. 

We denote the number of steps of the real-space length A in positive 
(negative) direction along the chain by ra^ (ra^) respectively. Obviously, 

nx = nl + n~, (67) 

and we also define 

dx = nl-n-, A = 1,2. (68) 
In the case considered, closed paths are those, for which 

c/i + 2^2 = 0. (69) 

To find r(ni,n2), let us first fix n^,n^ on some values permitted by (1691) . 
The number of closed paths with exactly n^,n^ 'positive' steps is 



riij \n{ J J n{ 



^The Reader is warned to not confuse the length of the 'step', A, with length of the 
path, n, being the total number of steps. 
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This is because we have to choose the order in which rii 'short' (A = 1), and 
n — rti 'long' (A = 2) steps appear in the path, and then choose those n]*" (n^) 
of the rii (^2) steps of each type, that are taken in the respective positive 
direction. To obtain the ri(ni,n2) we have to sum (1701) over allowed values 

of Til , . 

Before we do that, it is convenient to express first nj^,n^ through nx,dx 
with the help of the equations (!U7l EH]), and utilize the constraint to 
eliminate di. The only difficulty is that permitted values of d2 depend in a 
nontrivial way on ni,n2. To find the explicit form of this dependence, we 
have to consider two cases, ni > 2n2 and rii < 2n2- In the first case, for any 
value of (i2 (M2I < ^^2)5 we are able to find such a corresponding di that the 
condition flU^ holds. In the case rii < 2n2, value d^""^ depends on rii, and 
also on the parity, V, of both ^rii and ^2, in the following way 

^max ^ if p^!^^ _ V{n2), - 2) otherwise. 

Note, that for rii = 2n2 results for both cases coincide. Finally, we obtain 
the following formula for r{nin2), 

^ ^ riM ^\-2{n^-2d2))\\{n2-d2))' ^ ' 
and after simple algebra, for the coefficients of the expansion fl66l) . 

^ 5 (i(ni + 2d2))\{\{n, - 2rf2))!(|(n2 + t^2))!(i(n2 - ^2))! ' ^^^^ 
Summation in ( 1721) is over the elements of the set V2, 

V2 = {-^r^ -dT^ + 2, . . . , rfr"" - 2, c/r""}, 

and we recall once more the dependence of d™"^ on ni,n2, 

• for Hi < 2n2 

^max ^ if p(^r^) ^ p^^^)^ _ 2) othcrwlse, 

• for ni > 2n2 d^"-^ = n2. 

Equation flMj) . combined with the explicit form of I'ni,n2! fT^l) . is a central 
result of this section, and may be regarded as generalization of ( !63|) . 

The coefficients i^ni,n2 have a quite complicated structure, however, there 
exist some relations between them. This is because, due to identity 

cos(2A;) = 2cos^(A;) - 1, (73) 
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Zl{(,i,(,2) obeys the following partial differential equation 

This implies the relations between the coefficients (172 p . 

(rii + 2){ni + l)L„^+2,n2 - (^2 + l)L„^,n2+i - 2L„^,n2 = 0. (75) 

2.3 Triangular and bcc lattices 

In this section we examine two-dimensional triangular, as well as three- di- 
mensional bcc lattices. The reason why we consider those two cases together 
is that the formulas giving the number of closed paths have a very similar 
form for both lattices. This is due to the similar geometry of triangular and 
bcc lattices, with the nearest neighbors of a given site forming two interpen- 
etrating simplices in two (three) dimensions, respectively. 



2.3.1 Triangular lattice 

Each site of an infinite triangular lattice has six nearest neighbors grouped 
in pairs along three lattice directions, their positions are given by ±ej, 

ei = (l,0), e2 = i(-l,v^), e, = ~{l,V3). (76) 

The reciprocal lattice of a triangular lattice is also a triangular one. The first 
Brillouin zone is a hexagon centered in the origin, with the vertices 

±(27r,-^), ±(27r,+^), ± (47rv^/3, 0)), (77) 

and area equal sVSn'^. In the Hamiltonian ([71) we assume nonzero hopping 
integrals only between the nearest neighbors of a given site i. In general, we 
can relate a different hopping integral to each of the three lattice directions. 
Instead, we concentrate on the symmetric case ti = ^2 = ^3 = t, there is 
however no conceptual difficulty to extend the following results to the most 
general one. The dispersion relation reads now 



exikx^ky) = 2tcos{kx) +4tcos ^"^^ '^'^^ 
and, consequently. 



VSky 



(7J 



2cos(fc:r)+4cos(^)cos( '^'"^ 
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Now we pass to the problem of finding ^t(0 coordinate representation. 
We denote the number of steps in the direction given by ±ej as ±nj, with 
nf + n^ = Hi, and define 

di = nl-n^, i = 1,2,3. (80) 

Let A^. (Ay) be the total displacement along x (y) coordinate. We have then 



1 

Ar, = di - -{d2 + ds). Ay = — ((^2 - 4), 



and, as for a closed paths we require A^; = Ay = 0, we obtain the condition 

d^ = d2 = d3 = d. (82) 

One can easily convince himself that maximal allowed value of d is 

d^ax = min(?ii, na, rig). (83) 

The number of closed paths with exactly nf,n~ 'positive' and 'negative' 
steps in each direction is equal 



n\ ni\ n2\ ris 



and the total number of all paths of length n is then obtained by summing 
first over all permitted values of n^, and then over all values Ui = nf + 
obeying 

Til + n2 + ns = n. (85) 

It is, however, convenient to eliminate first nf{n~) with the help of ( l80l) . and 
apply the conditions ( l82l) . (l83i) . As the result, the coefficients T„ = T{n)/n\ 
are given by 



constr 



n\ ^ ni\n2\n^\ ^ \\{ni + d) J \\{n2 + d) J \\{n^ + d) J 

The first summation in ( l86i) is taken over those ni,n2 and which satisfy 
condition ( 185|) . What is more, ni,n2 and ng must have the same parity, for 
even n all are even, whereas for odd n all rii must be odd, too. This 
reflects the fact, that any closed path can be regarded as combination of 
a number of elementary paths of the two kinds: first, with n = 2, is a 
step forward and step backward in some direction, the second, with n = 3, 
consists of one step in each direction. The second summation runs over 
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d e {-dmax, -dmax + 2, . . . , d.max " 2, d.max} ■ In formula ([HnD some factorials 
cancel, and we finally obtain 

Explicitly computing ^t(0 np to sixth order using the above formula, we 
obtain 

zt{o = i+2e+ 2e + + + + o{e). (ss) 

4 o 

This coincides with the result obtained by integration of the expansion of the 
exponent in (1791) to the same order in ^. 

2.3.2 Bcc (body centered cubic) lattice 

This case is in many respects very similar to the just examined triangular 
lattice, with the main difference being higher dimensionality of the bcc lat- 
tice. Each lattice site has now eight, instead of six, nearest neighbors, given 
by vectors ±ej, 

ei = ^(1,1,1), 62 = ^(-1,-1,1), 

62 = ^(-1,1,-1), 64 = ^(1,-1,-1). (89) 

Similarly to the case of triangular lattice, non-zero hopping integrals are 
assumed only between nearest neighbors, and we put all of them equal, ti = 
h = h = = t (consequently, = 1^2 = ^3 = ^4 = 0- Reciprocal lattice 
of the bcc lattice is the fee lattice (see [3], p. 74), the first Brillouin zone 
is regular rhombic dodecahedron with the volume V = IGvr^. The dispersion 
relation is k k k 

GbcciK, ky, k^) = 8tcos(-^) cos(-^) cos(-^), (90) 

and the partition function in k-representation reads 

^bcciO = / e««™^('^)™^(^)™^(^)rfA;,rfA;,rfA;, = ^nC- (91) 

The coefficients Bn are obtained in a way very similar to those for the trian- 
gular lattice, also the notation is analogous and self-explanatory, we skip the 
derivation then and present only final results. Bn are given by 
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The first summation in fl92l) is taken over those rii, ?t,2, ^3, 77,4 which satisfy 
the condition rii + 77-2 + ^3 + 77,4 = n, where rij is a number of steps in the 
direction given by ±ej. Likewise in the previous section, all Ui must have 
the same parity, which implies that n must be even in the present case. The 
second summation runs over the elements of the set 

T-^bcc { dmaxi dmax ~l~ 2, . . . , d^aax 2, dmax} (^3) 

where 

dmax = min(?7i, n2, n^, 774). (94) 
After obvious simplifications, we finally obtain 

^ ^^ + dmiin, - d)y: ^^^^ 

Explicit evaluation of ZbcdO ^'^ twelfth order in ^, either using the formula 
( l95l) or by numerical integration of ( l9Ti) . gives the expansion 

e(o = 1 + 9e + + + '^^r + + an. (96) 

Interestingly, for n < 30 all i3„ are squares of rational numbers. The question 
arises, is this true for all values of n? However, at the moment we are not 
able to prove this conjecture or to find counterexample. 



2.4 Honeycomb and diamond lattices 

Likewise in the previous section, here we also analyze together two lat- 
tices, namely the two-dimensional honeycomb (graphene) lattice and three- 
dimensional diamond lattice. Once again, the reason is the similar geometry 
of both lattices, in the present case three (four) nearest-neighbors of a given 
site form a two (three) dimensional simplex, respectively. This results in a 
very similar form of the Taylor series expansion for the partition function in 
this two cases. Both honeycomb and diamond are non-Bravais lattices, with 
the basis consisting of M = 2 atoms, we denote the resulting two sublattices 
A and B, respectively. As a consequence, evaluation of the partition func- 
tion in quasimomentum representation is not that straightforward as in the 
Bravais (M = 1) case. However, this is a standard textbook problem, so we 
omit the detailed explanations. What is important, our method applies here 
without any serious modification. 



21 



2.4.1 Honeycomb lattice 



The nearest neighbors of any lattice site j G A of an infinite honeycomb 
lattice are given by vectors 

"'-(0.—). <=^=(i.^). <'3=(-^.f)- (W) 

whereas the neighbors of site i G -B by — ei, — e2, —e^. Because honeycomb 
lattice is, in fact, a triangular Bravais lattice with two-site atomic basis, its 
reciprocal lattice is also a triangular lattice, identical to that considered in 
section 5.1, due to the particular choice of vectors ( I971) . Consequently, also 
the first Brillouin zone is the same, i.e. it is a hexagon given by ( 1771) . 

I the Hamiltonian ([71) the non-zero values of the hopping integrals are 
assumed only to three nearest neighbors of a given site, and we put them 
equal. In the present case the dispersion relation consist of M = 2 subbands. 



\ 



3 + 2 cos(A;^) + 4 cos ( y cos I I , (98) 



where a = — 1 (a = 1) refers to lower (higher) subband, respectively. To eval- 
uate the partition function we have to integrate states within each subband, 
and than to add contributions from both of them. Then 



ZgIO = — 1= — / e V 2 \^ 2 y^L 



O-g. /3+2cOs(fci:)+4 COSi ^) COS 



TC^ J FBZ 



(T=-l,l 



= -4—/ cosh{^eG,i{h,ky))dhdky = J2GnC- (99) 

Now let us evaluate Zg{C,) in Wannier basis. We have S = MA = 2A, 
terminal vertex j of any closed path may belong either to sublattice A or 
B. Those two cases are equivalent with respect to the way we compute the 
number of closed paths, the only consequence is an additional factor 2. 

In order to close the path, the number of steps in each of three lattice 
directions, taken form A to B must be equal to that of steps form B to A, 
consequently a path must be of the form A^B— >-A— >■...— >-B^ A. 
Denote the number of steps, taken from j G A in the i-th direction (i.e. given 
by ±ej) as pi. To find r(n) we have to choose those pi of the total p steps 
A ^ B that are taken in i-th direction, and then independently do the same 
for steps B ^ A. The resulting number is equal (p!/(pi!p2'P3')) • Summing 
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over all values oi pi,p2,P3, which in the present case are of arbitrary parity, 
but obey the constraint pi + p2 + P3 = p = n/2, we obtain 



where the factor 2 was included according to the previous discussion. Using 
formula f llOOp we find the coefficients Gn for n < 6. This gives us 

G(o = 2 + + + + o{e). im 

This is precisely what one gets expanding ( l99l) up to sixth order in ^ and 
doing the remaining integrations. 

2.4.2 Diamond lattice 

Diamond lattice is related to the honeycomb lattice in very much the same 
way as bcc to the triangular one. It can be regarded as a fee (face cen- 
tered cubic) lattice with two-atom basis, and thus two sublattices, A and B. 
Consequently, the first Brillouin zone is that of fee lattice, namely truncated 
octahedron, its volume is V = 327r^ (see [3], p. 76, with a = 1). Each site 
of sublattice A [B) has four neighbors, their positions are given by (— e^), 
respectively, where ej vectors are the following 

61 = ^(1,1,1), e2 = ^(-1,-1,1), e2 = ^(-1,1,-1), e4 = ^(1,-1,-1). 



(102) 



Likewise in the honeycomb lattice case, the dispersion relation consist of two 
subbands labeled by cr = — 1, 1; it is related to the dispersion relation of the 
fee lattice as follows, 



eDaiK, ky, kz) = cTt\J 4 + Aefcc{K, ky, k^) (103) 

where 

€fcc{kx,ky,kz) = 4t^COS yCOSy + COS y COS y + COS y COS y (104) 

The resulting partition function is then 

= T7~^ I cosh(^A/4 + 4e/cc(fcx, ky, kz))dkxdkydkz. (105) 
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The way we compute the number of closed paths in the present case directly 
follows the case of honeycomb lattice, with the only difference that now we 
have three, instead of two, spatial dimensions. We have then 

Vn{n) = 2 J2 ( , f , , )' = ^!gn, (106) 

rii = 2pi, i = 1, . . . ,4, n = 2p, and sum over pi is constrained, J2t=i Pi ~ 
but again each pi may be of arbitrary parity, and the factor 2 appears due 
to the number of equivalent sublattices. ^d(0 expanded up to eighth order 
in ^, reads explicitly 

2 + + le + %e + + oie')- im 

3 Appendices 

3.1 Appendix A. One to one correspondence between 
the terms given by the formula (1301) and respective 
closed paths 

In this paragraph we give a simple justification of the fact, that there is a one 
to one correspondence between terms (1301) and closed paths they generate. 

Clearly, any term (IHUl) leads to some path. We have to show that different 
such terms correspond to different paths. Let us consider some term in 
question, 

A. . A A. . A- ■ A A- ■ A- ■ A A- ■ 

(108) 

and assume, that apart from the above, there exist some other ordering of 
its n Aij operators, which gives the same closed path. Any such ordering 
may be obtained by some permutation of Ai^ operators in formula (llOSp . it 
is enough then to examine the case of a transposition. Let us swap ^ij^i^.i 
and Ai.i._^, assuming first that Ai.i._^ 7^ Ai^i^ and Ai^i^^_^ ^ A^i^, 

A. . A A- ■ A- ■ A A- ■ A- ■ A A- ■ 

^lUn^lnln-l ^lk+1 Ik ^Ij Ij- 1 ^Ik- 1 lk-2 + 1 Ij ^Ifclk- 1 ^Ij - 1 lj-2 ^\1W 

(109) 

Obviously if A^ik-i — this transposition does not affect (11081) . yield- 

ing identical term. In the opposite case, after writing some A\^ operators in 
terms of creation and annihilation operators, (11091) reads 

^iiin^inin-i ' ' ' ^ik+1 ifc ^ij ^ij- 1 ^ik- 1 ik-2 ' ' ' '^ij+i'^ij^ikik-i^ij-iij-2 ' ' ' ^i2ii- 

(110) 
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There are two possibilities now: either some operators in between q. and 
Ci. carry site index ij or this index does not appear there. In the former 
case, after commuting some operators, a part of the term f llOQp has the form 
• • ■ cl,{l—ni.)ci. ■■■ = ■■■ c|.Cij ■•■ = ■■■ in the latter case that part also 

reduces to ■ ■ ■ cl,c\.- ■ ■ = ■ ■ - rii.- ■ ■ . Consequently, in both cases the whole 
term fll09p vanishes, due to the fact that we compute trace using one-particle 
states only. The same reasoning applies to the reversal of order in which C;^ 
and Ci^ appear in f llOQp . 

If one, but not both, of the permuted operators is A^^i^ or Ai^i^ the above 
reasoning applies without a modification. Finally, if we interchange Ai^i^ 
with in order to obtain the non-zero term, 12 = in must hold. But then 

we obtain the path with terminal vertex 12 7^ ii, so different from our point 
of view. 

3.2 Appendix B. Linear chain of atoms with complex 
nearest-neighbor hopping integrals 

Consider the case analyzed in section 12.11 but let the hopping integral be an 
arbitrary complex number, t = \t\e^'^. This leads to a modified dispersion 
relation, 

ec{k,ip) = 2\t\ cos {k + (p). (Ill) 
Consequently, the partition function is given by 

Zc{p,v,A) = j exp(-2pcos(A; + (^)), (112) 

keFBZ 

where ^ = —pt = — pe**^ and p = \^\. Zc{p,'^,A) is obviously real, periodic 
function of (p with a period 27r. Expanded in a Fourier series it reads 

00 

Zcip,if,A) = C,{p)e^P^, C;{p) = C.,{p). (113) 

p=— 00 

On the other hand, from the way we compute in coordination represen- 
tation, it is clear that each step in a positive (negative) direction brings a 
factor e**^ (e"**^), respectively, and that each path is characterized by a total 
phase ipd, 

This shows, that in the expansion f lllSp the only nonzero terms are those 
with p = d = cA, and d given by ( l55l 1581) . Comparing (11140 and (11130 . we 
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see that Cd{p) = C^d{p), and that (11131) can be rewritten as 

^£(P, V^, A) = ao(p) + ^ ad{p)cos{(pd), (115) 



d{n)>0 

where 



keFBZ n=0 {\2"')-) 



k&FBZ n=0 1-2 V /J L2V 

(117) 

Assume now even A, then n = 2z/, d = 26. Equating flll4p to (11151) and 
having in mind (I116[ I117p we obtain the following generalization of ( I6T1) . 



Eexp(-2pcos(/c + (y9)) / 1 9 cos(25v9) \ 

A " ^\^.^ ^ \(u + 6)V.\u-6]\) 

keFBZ u=0 <5(i/)>0 '-^ ^' ^ ' 



u=0 S{u) ' /J L J 



If we put Lf = TT in the above formula we obtain the results given by 
On the other hand, for (y9 = | we have cos{26ip) = cos(57r) = (—1)*^ and then 

- Y: exp(2ps.nW)=^(^ ' )p'-. (119) 

keFBZ u=0 5(iy} '-^ ^' ^ ' 



Summary and conclusions 

In this paper we have presented the method of evaluating the partition func- 
tion of a single-electron, periodic system, using specific basis, in which the 
system Hamiltonian is not diagonal. Utilizing the symmetry properties of 
the system and its Hamiltonian, as well as invariance of trace operation with 
respect to change of the basis in the Hilbert space, we are able to establish 
several mathematical identities. They could be found important from the 
point of view of both pure and computational mathematics. 

The possible extensions of the present works include first analysis of infi- 
nite lattices not examined yet, like fee or kagome lattice, as well as detailed 
investigation of finite lattices for geometries more complicated than that of 
linear chain, and also for various PBC. The extension of the method to many- 
electron systems or systems without translational invariance is also possible 
in principle, but the technical difficulties that arise make the application of 
combinatorial techniques very problematic in these cases. 
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